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The magnetohydrodynamical behavior (MHD) of accretion disks is reviewed. 
A detailed presentation of the fundamental MHD equations appropriate for 
protostellar disks is given. The combination of a weak (subthermal) magnetic 
field and Keplerian rotation is unstable to the magnetorotational instability 
(MRI), if the degree of ionization in the disk is sufficiently high. The MRI 
produces enhanced angular momentum and leads to a breakdown of laminar 
flow into turbulence. If the turbulent energy is dissipated locally, standard 
"a" modeling should give a reasonable estimate for the disk structure. Be- 
cause away from the central star the ionization fraction of protostellar disks 
is small, they are generally not in the regime of near perfect conductivity. 
Nonideal MHD effects are important. Of these, Ohmic dissipation and Hall 
electromotive forces are the most important. The presence of dust is also crit- 
ical, as small interstellar scale grains absorb free charges that are needed for 
good magnetic coupling. On scales of AU's there may be a region near the 
disk midplane that is magnetically decoupled, a so-called dead zone. But the 
growth and settling of the grains as time evolves reduces their efficiency to 
absorb charge. With ionization provided by coronal X-rays from the central 
star (and possibly also cosmic rays), protostellar disks may be sufficiently 
magnetized throughout most of their lives to be MRI active, especially away 
from the disk midplane. 



1 Introduction 

An understanding of the magnetohydrodynamics (MHD) of protostellar disks 
is crucial for the theoretical development of these objects. There is no get- 
ting around the fact that the subject is decidedly inelegant. In principle, we 
are only interested in solving the equation F = ma, but the range of top- 
ics that bears on this problem is truly daunting. Molecular chemistry, dust 
grain physico-chemistry, photo-ionization physics, aerosol theory, and non- 
ideal MHD are all key players in this game. If the results of the theorists' 
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efforts had been little or no progress, there would have been no shortage 
of excuses. ("More work is needed.") But in fact there has been substantial 
progress in understanding important issues, more perhaps than some practi- 
tioners may even realize. A sort of consensus is beginning to take shape of 
the gross properties of protostcUar disks that arc dictated by the demands of 
MHD physics. In this review, I will try to put forward as strong a case as I can 
for what I somewhat boldly regard as the canonical protostellar disk model, 
while at the same time try not to gloss over what are genuine uncertainties 
or difficulties. 

Protostellar disks are a class of accretion disks, and one area where there 
certainly has been a great deal of progress in recent years is in the develop- 
ment of accretion disk theory. The realization that a combination of differen- 
tial rotation and a weak magnetic field is profoundly unstable and produces 
turbulence has given the subject a foundation on which to build (Balbus & 
Hawley 1991, 1998). Indeed, the primary reason for studying MHD processes 
in protostellar disks in detail is that once these systems are no longer self- 
gravitating, magnetic fields profoundly influence their dynamical behavior. 
It is not possible to understand angular momentum transport in low mass 
("T-tauri") disks without focusing on magnetic fields. 

Magnetic fields provide a conduit for free energy to flow from the differen- 
tial rotation to the disk itself, producing fluctuations, turbulent heating, and 
quite possibly the directly observed outflows. Most importantly, we shall see 
that MHD turbulence produces well-correlated fluctuations of the radial and 
azimuthal components of the magnetic field and fluctuating gas velocity. It is 
precisely these correlations that result in a substantial outward transport of 
angular momentum, allowing the disk gas to spiral inward and accrete onto 
the central star. 

It may also be possible in principle for winds to remove angular momentum 

from the disk without completely disrupting it. (Sec Kocnigl, this volume.) 
The difficulty that needs to be overcome for this mechanism is that the out- 
flow must occur over the whole face of the disk, not just from the innermost 
regions from where jets are launched: material near the inner disk edge has 
already lost almost all of its angular momentum. Instead, the outer disk ma- 
terial, torqued by the magnetic field, must continuously slip relative to the 
field, or else the field tends to become very centrally concentrated (Lubow, 
Papaloizou, & Pringle 1994). How this might work without generating inter- 
nal disk turbulence that itself transports angular momentum remains to be 
fully understood. 

This understanding that magnetic fields play a crucial role as the source 
of disk turbulence developed more than thirty years after the basic instability 
(now known as the magnetorotational instability or MRI) first appeared in 
the literature (Velikhov 1959). But even if the importance of the instability 
had been immediately grasped, without the computational power that became 
available only post 1990, it is doubtful that the impact would have been the 
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same. There is no substitute for being able to visualize on your desktop the 
development of a linear instability into full-fledged MHD turbulence. 

The MRI depends upon the presence of electrical currents to do its job, 
and this means that the disk gas must be at least partially ionized. "Partially 
ionized" in practice could mean even a minute electron fraction. Tiny traces of 
electrons can magnetize the disk, an eflfect we will quantify in §8. It is because 
even a wisp of a magnetic field and the merest trace of electrons take on such 
dynamical significance that protostellar disk dynamics depends so heavily on 
protostellar disk chemistry. 

Young protostellar disks need sodium, potassium and other trace "vita- 
mins" to make them vigorous and active, because these alkalis arc important 
in regulating the delicate ionization balance in the gas. It is ultimately the 
ionization level that determines the resistivity that determines whether the 
magnetic field is well-coupled to the gas or not. This is not an issue in accretion 
disks around compact objects, for which even a small fraction of the turbu- 
lent energy dissipated would be enough to ensure that the gas is thermally 
ionized. Protostellar disks, on the other hand, are big (so the free energy of 
differential rotation is small), dusty (so grains absorb free electrons), and cold 
(so thermal ionization is unimportant except near the star). The overarching 
question of protostellar MHD research is to understand under what conditions 
the abundance of free electrons level drops below the level needed to ensure 
good magnetic coupling. The magnetic coupling leads in turn to instability, 
turbulence, and enhanced transport, the essence of disk dynamics. Much of 
this review will be devoted directly to the question of disk ionization balance. 
For readers wishing additional astronomical background material, the text of 
Lee Hartmann (1998) is an excellent choice. Balbus & Hawlcy (2000) and 
Stone et al. (2000) are earlier reviews of MHD transport processes in proto- 
stellar disks. More general disk reviews include Pringle (1981; a classic review 
but pre- MRI), Papaloizou & Lin (1995), Lin & Papaloizou (1996), Balbus & 
Hawley (1998), and Balbus (2003). In addition, the Protostars and Planets 
series published by the University of Arizona Press series contains a useful 
historical record of the development of the subject of protostellar disks; the 
latest volume is Protostars and Planets V (Reipurth, Jewitt, & Keil 2007). 

2 On the Need for MHD 

The modern view of protostellar disks is heavily influenced by accretion disk 
formalism. Lin & Papaloizou (1980) is the first pioneering study to bring to 
bear accretion disk formalism to the study of protostellar disks, though it is 
now thought that MHD turbulence, rather than convectivc turbulence, is the 
key to protostellar disk dynamics. To understand the central role of a magnetic 
field in our understanding of protostellar disks, it is of interest to review where 
we stand with respect to our knowledge of the onset of turbulence in rotating 
flows. 
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The principle problem for classical disk theory (Shakura & Sunyaev 1973) 
was to discover the origin of the very large stresses that were needed to trans- 
port angular momentum, a process then (and still often now) referred to as 
"anomalous viscosity." The approach of Shakura & Sunyaev was to make a 
virtue of necessity, arguing that the small microscopic viscosity implied a very 
large Reynolds number (Re)0, and a large Reynolds number meant that shear 
turbulence was present. This conclusion was sustained by the belief that the 
laboratory experiments showing a nonlinear breakdown of Cartesian shear 
layers at values of Re in excess of ~ 10^ would naturally carry over to Keple- 
rian disks, where its value was considerably higher. The fact that the inertial 
force associated with the Coriolis effect (not present in planar flow) was well 
in excess of the destabilizing shear does not seem to have been viewed as 
troublesome to proponents of this view. 

The classical laboratory method for investigating the stability of rotating 
flows is a Couette apparatus. In such a device, the space between two coaxial 
cylinders is filled with a liquid, almost always ordinary water. The two cylin- 
ders rotate at different angular velocities, let us say i7i„ and ilout- The gap 
between the cylinders is of order a centimeter, and a stable rotation profile 
will be attained in a matter of minutes by viscous diffusion, even if Re is very 
large. (An unstable rotation profile will of course never be found as such; the 
flow will remain permanently disturbed.) By choosing J7i„ and /^out appro- 
priately, a section of a Keplerian disk can be accurately mimicked, and the 
effects of Coriolis forces can be studied. 

According to the classic text of Landau & Lifschitz (1959), a Couette ro- 
tation profile that is found in the laboratory to be stable "does not actually 
mean, however, that the flow actually remains steady no matter how large 
[Re] becomes." The monograph by Zel'dovich, Ruzmaikin, & Sokoloff (1983) 
is even more explicit on the question, unambiguously stating that laboratory 
experiments had already shown that Keplerian rotation profiles were nonlin- 
early unstable at sufficiently large Re. 

These are stunning claims. At the time these statements were written, 
there was certainly no credible laboratory evidence that Keplerian rotation 
profiles were nonlinearly unstable. The first explicit claim that Keplerian flow 
was unstable based on laboratory evidence seems to be the unpublished result 
of Richard (2001). However, a later experiment by Ji et al. (2006) suggests that 
the earlier finding of instability was probably due to spurious effects arising 
from the interaction between the flow and the endcaps of the experiments 
(Ekman layersH). When care is taken to minimize this interaction, Keplerian 
flow is found to be linearly and nonlinearly stable at values of Re up to 

^ The Reynolds number Re is the ratio of a characteristic flow velocity V times a 
charateristic length L divided by the kinematic viscosity, u, Re = VL/v. 
Ekman layers are narrow fluid layers between a solid boundary and a rotating 
flow. When the bulk of the flow is not at rest in the frame of the boundary, the 
fluid rotation changes rapidly in the Ekman layer and can produce a secondary 
circulation pattern invading the entire flow. 



Magnetohydrodynamics of Protostellar Disks 5 



2 X 10^. Serious students of disk theory would do well to develop a deep 

skepticism of "large Re means disks arc turbulent" arguments, which sadly 
are still promulgated in contemporary textbooks. 

Why is the Coriolis force is so harmful to the maintenance of turbulence? 
Moffatt, Kida, & Ohkitani (1994) have referred to stretched vortices as the 
"sinews of turbulence." Turbulence is maintained in shear flows by vortices 
that are ensnared along the axis of strain, coupling the free energy of the shear 
directly to the internal vortex motion. Two neighboring fluid elements in a 
vortex are rapidly pulled apart as its circulation rises. This is possible only if 
Coriolis forces are absent. Despite the presence of shear, Coriolis forces induce 
epicyclic, oscillatory motions that do not allow anything resembling continu- 
ous vortex stretching. Without this, turbulence cannot be maintained. This 
argument breaks down of the strain rate much exceeds the formal epicyclic 
oscillation period, even the latter is quite well defined. Such flows, in fact, arc 
found to be nonlinearly unstable in numerical simulations. But for Keplerian 
flow, the epicyclic frequency is just the local O of the disk, and it exerts a 
strongly stabilizing influence. 

Figure (1) illustrates this point. Locally, Cartesian shear flow and Keple- 
rian differential rotation appear to be similar. But perturbed fluid elements 
behave very differently in the two systems. Viewed from a frame moving with 
the undisturbed flow, a displaced fluid element in shear flow approximately 
follows an unbounded parabolic trajectory A, while a displaced fluid element 
in a Keplerian disk approximately follows a bound epicycle, B. A small per- 
turbing vortex would be stretched continuously in shear flow, tapping into the 
free energy source needed to maintain turbulence. By contrast, the embedded 
vortex would merely oscillate in a disk. It is this difference, due entirely to 
the presence of the rotational Coriolis force, that is ultimately responsible for 
the hydrodynamical stability of Keplerian disks, whereas shear flow can be 
nonlinearly imstablc. 

The title of the Ji et al. (2006) paper is "Hydrodynamic turbulence can- 
not transport angular momentum effectively in astrophysical disks," and that 
seems as good a onc-scntence summary on the topic as any; their experi- 
ment is probably the final word on nonlinear Keplerian hydrodynamical shear 
instability. It is gratifying, therefore, that all of the difficulties encountered 
by Coriolis stabilization vanish when magnetic fields are taken into account. 
This is due to the fact that the magnetic field introduces new modes of re- 
sponse (shear waves) that are much less prone to Coriolis stabilization. One 
of these waves, the so-called "slow mode," turns into a local instability when 
differential rotation is present with a weak magnetic field. This is the magne- 
torotational instability. To understand this porcess in more detail, it is best 
to begin with the fundamentals of MHD. 
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Fig. 1. The stabilization mechanism in a rotating disk versus Cartesian shear flow. 
In shear flow, displacements are unbounded, as indicated by path A. In a disk, 
epicyclic motions keeps the displacement tightly bound. See text for further details. 



3 Fundamentals 

In this section, a detailed derivation of the fundamental MHD equations is 
presented. The discussion will be more technical here than in most of the rest 
of the paper, but it is very important to see how the basic governing equations 
of the subject arise, and much of this material is not so easy to find outside 
of specialized treatments. I hope the reader will have the patience to read 
carefully through this section. 

A protostellar disk is a gas of neutral particles (predominantly H2 molecules), 
electrons, ions (the most important of which will generally be K+ and Na'^), 
and dust grains. We defer to §8 a discussion of the complications due to the 
presence of dust grains and consider here the gas dynamical equations for 
mixture of neutrals, ions and atoms. 

Each species (denoted by subscript s) is separately conserved, and obeys 
the mass conservation equation 



^ + V.(p,t;,)=0 



(1) 



Magnetohydrodynamics of Protostellar Disks 7 



where Ps is the mass density for species s and Vg is the velocity. The flow quan- 
tities for the dominant neutral species will henceforth be presented without 
subscripts. 

The dynamical equation for the neutral particles is 
dv 

p—+p{vV)v = -VP-pV^-Pr,I-Pr,^ (2) 

ot 

where P is the pressure of the neutrals, <P the gravitational potential and 
Pni (Pne) is the momentum exchange rate between the neutrals and the ions 
(electrons). 

Let us examine these last two important terms in more detail. Pni takes 
the form 

Pni = np,.ni{v - Vi)v.rij (3) 

where n is the number density of neutrals, and /i„/ is the reduced mass of an 
ion and neutral particle, 

m/m„ 

A*n/ = \ , (4) 

TO/ + TO„ 

TO/ and TO„ being the ion and neutral mass respectively. Vni is the collision 
frequency of a neutral with a population of ions, 

Vnl = ni{anlWnl)- (5) 

In equation ([5]), nj is the number density of ions, cr„/ is the effective cross 
section for neutral-ion collisions, and Wni is the relative velocity between a 
neutral particle and an ion. The angle brackets represent an average over a 
Maxwellian distribution function for the relative velocity. (The mass appearing 
in this Maxwellian will of course be the reduced mass /x„/.) For neutral- ion 
scattering, we may take the cross section (t„/ to be approximately geometrical, 

— 1/2 

which means that the quantity in angle brackets will be proportional to /i^j . 
The order of the subscripts has no particular significance for the cross section 
(T„/ reduced mass /i„/ or relative velocity Wni- But i^jn differs from lyni- the 
former is proportional to the neutral density n, not the ion density nj. 
Putting these definitions together gives 

Pnl = nniPni{<TnlWnl){v - Vi) (6) 

In accordance with Newton's third law, this is symmetric with respect to 
the interchange n I, except for a change in sign, Pni = pin- All of 
these considerations hold, of course, for electron-neutral scattering as well. 
Explicitly, we have 

Pne = nneP.„e{o-neWne)iv - V^) ~ rWeTOe (cr„eti;„e) (■" " t'e)- (7) 



The gas will be locally neutral, so that rig = Zrii where Z is the number of 
ionizations per ion particle. In a weakly ionized gas, Z = 1. The reduced mass 
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Une may safely be set equal to the electron mass mg. We will use the following 
expressions for the collision rates (Draine, Roberge, & Dalgarno 1983 jl: 

{cTnlWnl} = 1-9 X 10"^ cm^ s^^ (8) 

10"^''^ (128fcT/97rTOe)^^^ = 8.3 x iQ'^^t^/^ cm^ s^^ (9) 



^ ne ^ne/ 



The electron-neutral collision rate is essentially the ion geometric cross sec- 
tion times an electron thermal velocity. (The peculiar factor of (128/97r)^/^ 
comes from the way the electron velocity is projected, when averaged over the 
electron Maxwellian distribution function, along the direction of the mean 
ion flow.) But the ion- neutral collision rate is temperature independent, much 
more beholden to long range induced dipole interactions, and significantly 
enhanced relative to a geometrical cross section assumption. Even if the ion- 
neutral rate were determined only by a geometrical cross section, |pti/| would 
excee \pne\ by a factor of order {nie/ UniY^"^ ■ In fact, the dipole enhancement 
of the ion- neutral cross section makes this factor yet larger. 

In the astrophysical literature, it is common to write the ion-ncutral mo- 
mentum coupling in the form 

Pin = ppnivi - V), (10) 

where 7 is the so-called drag coefficient, 

mj + TO„ 

and we will use this notation from here on. Numerically, 7 = 3 x 10^'^ cm"^ 
for astrophysical mixtures (Draine, Roberge, & Dalgarno 1983). 
The dynamical equations for the ions and electrons are 

and 



+ Pivi-Vvi = - VP/ - piV<^ + Zeni [e +^xB^ Pin (12) 



Pe-g^ + PeVe-VVe = -VPe - PeV^ - eUe + X B j - Pen, (13) 

respectively. Throughout this paper, e will denote the positive charge of a 
proton, 4.803 x 10~^° esu. For a weakly ionized gas, the Lorentz force and 
coUisional terms dominate in each of the latter two equations. Comparison 
of the magnetic and inertial forces, for example, shows that the latter are 
smaller than the former by the ratio of the proton or electron gyroperiod to a 
macroscopic flow crossing time. Thus, to an excellent degree of approximation, 



Zeni 



— xB)-p7„-0, (14) 



^ A recent calculation of H-H"*" scattering by Glassgold, Krstic, & Schultz (2005) 
may imply a slightly higher value for the neutral-ion collision rate. 
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and 

-ene(£;+^xB) -per. = 0. (15) 
The sum of these two equations gives 

— XB=pin+Pen (16) 
C 

where charge neutrahty = Znj has been used, and we have introduced the 
current density 

J EE ene{vi - Ve). (17) 
The equation for the neutrals becomes 

dv J 
p— + pv-Vv = -VP - pV<P+ - xB (18) 
ot c 

Due to colhsional couphng, the neutrals are subject to the magnetic Lorentz 
force just as though they were a gas of charged particles. It is not the magnetic 
force per se that changes in a neutral gas. As well shall presently see, it is the 
inductive properties of the gas. 

Let us return to the force balance equations for the electrons: 



en, 



(£;+^xb) -per. = 0. (19) 



After division by — erie, this may be expanded to 

1 Til/ h' 

E+-[v + iv,- vi) + {vi -v)]xB + [{v, - vi) + {vi - v)] = 0, 

c e 

(20) 

where we have introduced the collision frequency of an electron in a population 
of neutrals: 

Ven = n{aneWne). (21) 

We have written the electron velocity in terms of the dominant neutral 
velocity v and the key physical velocity differences of our problem. It has 
already been noted that in equation (fT6|) . Pen is small compared with pin, 
provided that the velocity difference |t;e — f | is not much larger than ~ ?;|. 
In fact, when used in (|20p . the p^n term in equation (|16p may quite generally 
be dropped: 

-xB-pir,. (22) 
c 

It may also be shown that the final term in equation (j20p 

(vi - v), 

e 

which would then be proportional to J X B, becomes small compared with 
the third term 
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-K - Vi)xB, 

c 

also proportional to JxB. Both of these rather technical claims are justi- 
fied in detail in Appendix A. (In both cases, the incurred error is of order 
{me/ pinY^"^ ■) These simplifications allow us to write the electron force bal- 
ance equation as 

V ^ J X B J (J X B) X B 

E + -XB + - = 0, (23) 

C erigC (Tcond c^lPPi 

where the electrical conductivity has been defined as 

2 

C^cond = (24) 

The associated resistivity rj is (e.g. Jackson 1975) 

^=JZ^^ (25) 

which has units of cm^ s^^. Numerically (e.g. Blaes & Balbus 1994, Balbus & 
Terquem 2001): 

' n 

Up 



234 ( — ) cm2 s-i (26) 



Equation (|23p is a general form of Ohm's law for a moving, multiple fluid 
system. 

Next, we make use of two of Maxwell's equations. The first is Faraday's 
induction law: 



We substitute E from equation ([23|l to obtain an equation for the self- 
induction of the magnetized fluid, 



dB „ 



„ JxB (J X B) X B cJ 

vxB h -^^ 



eue c^ppi 



(28) 



It remains to relate the current density J to the magnetic field B. This is 
accomplished by the second Maxwell equation, 

ilj.VxB + if (28) 

c c at 

The final term in the above is the displacement current, and it may be ignored. 
Indeed, since we have not, and will not, use the "Gauss's Law" equation 

V-£; = 47re(Zn/ - ne), (30) 

we must not include the displacement current. In Appendix B, we show that 
departures from charge neutrality in V-£J and the displacement current are 
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both small terms that contribute at the same order: jt? . These must both be 
self-consistently neglected in nonrelativisitc MHD. (The final Maxwell equa- 
tion V'B = adds nothing new. It is automatically satisfied by equation 
(|?7|) . as long as the initial magnetic field satifies this divergence free condi- 
tion.) These considerations imply 

J=— VxB (31) 

4-77 

for use in equation (|28p . 

To summarize, the fundamental equations of a weakly ionized fluid are 
mass conservation of the dominant neutrals (eq. [T] ) 

|^ + V.(p^)=0, (32) 
the equation of motion (eq. [18] with [31]) 

+ /5i;-Vi; = -VP-pV^ + — (VxJ5)xB, (33) 
at Att 

and the induction equation (eq. [55] with [55] and [21]) 



dB 



„ c{V X B)xB UV X B) xB]x B _ „ 
vxB- -i '- + ^ i rjV X B 



(34) 

It is only natural that the reader should be a little taken aback by the sight 
of equation ((34| . Be assured that it is rarely, if ever, needed in full generality: 
almost always one or more terms on the right side of the equation may be 
discarded. When only the induction term vxB is important, we refer to this 
regime as ideal MHD. The three remaining terms on the right are the nonidcal 
MHD terms. 

To get a better feel for the relative importance of the nonidcal MHD terms 
in equation (j28|) , we denote the terms on the right side of the equation, moving 
left to right, as / (induction), H (Hall), A (ambipolar diffusion), and O (Ohmic 
resistivity). Then, the scaling of each of these terms relative to the Hall term 
i? is 

I V A UJcI O Ven 

— , -TT^^ ' 77 35) 

H VI - Ve H 7(0 H LOce 

where the ion cyclotron frequency ojd — eB/mjc, and similarly for the elec- 
tron frequency ujce, with me replacing m/. 

We will always be in a regime in which the presence of the induction term 
is not in question, i.e. the relative ion-electron drift velocity vj — Ve is always 
comparable to or in fact much less than v. More interesting is the relative 
importance of the nonideal terms. The explicit dependence of A/H and 0/H 
in terms of the fluid properties of a cosmic gas is given in Balbus & Terquem 
(2001): 
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cs 



(37) 



and 

H V8 X 10" cni-3 J \VA 

Here n is the total number density of all particles, T is the kinetic temperature, 
va is the so-called Alfven velocity, 

R 

(38) 

4:71 yL/ 

and C5 is the isothermal speed of sound, 

kT 

4 = 0.429 (39) 

Ulp 

where k is the Boltzmann constant and rUp the mass of the proton. The coef- 
ficient 0.429 corresponds to a mean mass per particle of 2.33mp, appropriate 
to a molecular cosmic abundance gas. 

As reassurance that the fully general nonideal MHD induction equation 
is not needed for our purposes, note that equations ([36]) and p7| imply that 
for all three nonideal MHD terms to be comparable, T 10^ K! Obviously 
this is not an issue for protostellar disks. In figure (2), we plot the domains of 
relative dominance of the nonideal MHD terms in the nT plane. Protostellar 
disks are dominated by the resistive and Hall nonideal MHD terms, except in 
the innermost regions, where ohmic dissipation is largest, and the outermost 
regions, where ambipolar diffusion becomes important (Wardle 1999). 

Our emphasis of the relative ordering of the nonideal terms in the induc- 
tion equation should not obscure the fact that ideal MHD is often an excellent 
approximation, even when the ionization fraction is ^ 1. For example, the ra- 
tio of the ideal inductive term to the ohmic loss term is given by the Lundquist 
number 

L = ^ (40) 
V 

where L is a characteristic gradient length scale. To orient ourselves, let us 
set L — O.IR, where R is the radial location in the disk. Then L is given by 

L ~ 2.5{ne/n){vA/cs)Rcm, 

Rem being the radius in centimeters. In other words, the critical ionization 
fraction at which L = 1 is about 

{ne/n)c„t = OA{cs/vA)R-r^ - 10-^^{cs/10va) 



at R — 1 AU. The actual nebular ionization fraction at this location may be 
above or below this during the course of the solar systems evolution, but the 
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point worth noting here is that Rem is a large number for a protostellar disk, 
whatever the ionization source! Ionization fractions far, far below unity can 
render an astrophysical gas a near perfect electrical conductor. It therefore 
makes a great deal of sense to begin by examining the behavior of an ideal 
MHD fluid. 



4 Ideal MHD 

The fundamental equations of ideal (single-fluid) MHD are mass conservation 

o 

^ + V.ipv)^0 (41) 
the dynamical equation of motion 

- + (..V)„.--v(p+-)-Vi.+ (-.)VB (42) 

where (p is the (central) gravitational potential function, and our newly sim- 
plified induction equation 

— = VxivxB) (43) 



14 Steven A. Balbus 



We shall work in a standard cylindrical coordinate system {R, (j), z), where R 
is the radius, (f> the azimuthal angle, and z the vertical coordinate. In these 
coordinates, the three components of the equation of motion are are 



— + wV 



1 d 



R 



pdR 



B 



-^^—Tr7,[P+^\-^ + ~VBR- 



Stt 



Im 



B 



Bl 



d „■ 



p + — 



B^ 
Sir 



Vz = - 



ld_ 
p dz 



P + 



B^\ d$ 



B 



Stt J 



dz 47rp ^ ^ 



and the three components of the induction equation are 

Br = -BrV-V + B-VVR 



d_ 

dt 



B. 



^ = V 

R R 



BV 



(I) 



dt^^-^ 



B, = -B^V-v + B-Vvz 



(44) 

(45) 
(46) 

(47) 

(48) 
(49) 



For many problems of interest, an important simplification can be made to 
these equations. The essence of rotational dynamics is local. Imagine going to 
a location Rq in the disk. The Keplerian angular velocity i? at i? = i?o is i?o . 
We hold J7o fixed but allow Rq to become arbitrarily large. Thus, V(j,{Ro) is 
unbounded, but v^/Rq = f2o is finite. 

Next, we introduce local coordinates {x, a, z) defined by 

R = Rq + x, x^Ro, (50) 

a = ^- Qot<^7r (51) 

and z is unchanged. We will formally introduce a new coordinate t' = t. This 
is desirable because we wish to take a time derivative at constant x, y, z, not 
at constant R, (j), z. In fact, 



dt'~dt^ 



(52) 



Partial derivatives with respect to R and (j) become partial derivatives with 
respect to x and a respectively. The Lagrangian derivative 



becomes 



9 „ 
dt+^-^ 



— + (t) - R^oe^yV 



(53) 



(54) 
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(The V operator now formally involves x and a derivatives.) This transfor- 
mation suggests that we introduce a new azimuthal velocity, 



(55) 



and for notational consistency we will use Wx and Wz for the radial and vertical 
velocities, though they are identical to vr and v^. 

In the local approximation, we assume that the magnitude w is much 
smaller than the (formally infinite) rotation velocity, and that w, the Alfven 
speed va, and the thermal velocity [P/ p)^/'^ are comparable in magnitude. In 
the limit Rn ^ oo, the radial equation becomes 



d_ 
dt' 



w-V 



wr - 2f2ow^ = R{f2l - Q\) 



ld_ 
p dx 



where the Keplerian angular velocity at R satisfies 



RO^ ^'^ 



Expanding to first order in R, 



'~dR 



The final form of the equation is 
w-V 



d_ 

di 



1 d 

Wx — 2Qwy = 3f2^x — 

pox 



P 



B2 



inp 

(56) 

(57) 

(58) 



B 

Anp 



•VBx, (59) 



where we have dropped the subscript and the prime '. The vector compo- 
nents X, y, z refer to the quasi-Cartesian system dx = dR, dy = Rda, and dz 
as before. (Sec fig. 2.) When non-Keplerian disks are considered, 3/?^ should 
be replaced by dfl"^ /dlnR. 

The transformed azimuthal {y) equation is 



9 „ 



1 d 



w,. 



and the z equation is 
'd 



dt 



+ w-'V 



dz 



pdy 



ld_ 
p dz 



2nwx = — + +1 — ^^v^ (60) 



Stt / 47rp 



P + 



52 \ 



Stt ) 



(61) 



The mass conservation and induction equations, like the z equation of motion, 
show no change of form in these corotating coordinates: 



1^ + V.{pw) = 0, 



(62) 
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— ^Vx{wxB) (63) 

To summarize, we have transformed our general equations into a coordinate 
frame rotating at the Keplerian angular velocity Qk{Ro), restricting the spa- 
tial domain to a small neighborhood around a particular fluid element. The 
rotational dynamics appears in the form of a Coriolis force —2flxw, and a 
centrifugal force that cancels the main gravitational force, leaving the residual 
tidal term iQ^x. These simplified local equations, which are the magnetized 
version of what is known as the Hill system, retain a rich dynamical content. 
This includes the full development of MHD turbulence. 



4.1 Magnetorotational Instability 

Let us apply equations ([56]) - ((63)) to the study of the paths of fluid elements 
departing from the origin x = y = z = Q. These Lagrangian displacements 
are in the xy orbital plane, and depend only upon z. We will ignore vertical 
stratification here, so that the equilibrium is z independent. Thus, we may 
assume a spatial dependence of e*'''^ for the perturbed fluid elements. 

We assume a very simple equilibrium magnetic field: a constant vertical 
component B, no components in the orbital xy plane. If the displacement 
vector of an element of fluid in the xy plane is ^, then the perturbed magnetic 
field at given location is 

(5J3 = Vx(|xB) (64) 

This result, which is true quite generally, is lengthy to prove by a direct 
assault, but can be intuited rather easily. The left side of the exact equation 
([63]) may be written 5B/6t, the change in B at fixed location divided by the 
change in time. In the equilibrium state this is zero. Imagine perturbing this 
pure state by giving each fluid an additional finite velocity U, but letting it 
act only for an infinitesimal time St. Then 

5B = Vx{w5txB) + Vx{U6txB) = Q + Vx{^xB) (65) 

where ^ = Udt is the displacement of the fluid relative to its equilibrium 
path. This result is perhaps more clear in integral form. If we integrate over 
a surface A and use Stokes Theorem, there obtains 

J 5B-dA = j {ix B)-ds = j {dsxi)-B, (66) 

where ds is a vector line element on the curve bounding A. This states that 
whatever explicit change in magnetic flux there is through A as the surface is 
displaced by ^, it is precisely compensated by the gain (or loss) of magnetic 
flux through of the side of the cylinder (elemental area ds X ^) swept out by 
A. This is just a round about way of saying that the magnetic flux through 
A is conserved when moving with the fluid itself. 
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For the problem at hand, 

6B = Vxi^xB) = ikxi^xB) ^ i{k ■ B)^, (67) 
since fc • ^ = 0. The magnetic term in equations ([56|l and ([60|) is 



We have introduced the so-called Alfven velocity 



VA = -5= (69) 



whose physical interpretation we discuss below. For the moment, note that 
the magnetic force is like tension or a spring: it is always restoring, and pro- 
portional to the displacement. 

In equations ([5S|) and (pO)) . the derivative 

D d „ 

is just the time derivative following a fluid element, so that we will write 
Dw/Dt = ^, and of course w = ^. Finally, we will allow for the possibil- 
ity of non-Keplerian rotation. The equations of motion for the fluid element 
displacements and i^j^ are then 



(70) 



iy + 2f2ix = -ikvAf^y (71) 

These are a simple set of coupled linear equations with time-dependent 
solutions of the form e'"'. Without loss of generality, we take w > 0. In the 
absence of rotation, the equations decouple completely and one finds 

iu = \kvA\ (72) 

These are waves that travel along the magnetic field lines with group and 
phase velocity va, exactly like waves on a string. These so-called Alfven waves 
are thus simple, nondipersive, incompressible disturbances. 
On the other hand, in the absence of rotation, we find 

is known as the epicyclic frequency. In most astrophysical applications, 
/2 decreases radially outwards and k < AfP. The fluid dispacements for this 
mode execute ellptical with a major to minor axis ratio of njlfl. The ma- 
jor axis lies along the azimuthal direction. The sense of the rotation around 
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the ellipse is retrograde, opposite to the sense of f2. These paths are known 
as epicycles. Epicyclic motion corresponds to the first order departure from 
circular orbits for simple planetary motion. 

The full dispersion relation resulting from equations (|70p and (j7ip is more 
complicated than one might have expected: 



d\nR 







(74) 



This is clearly not just a matter of adding Alfvenic and epicylic motion in 
quadrature, something new is going on. The dispersion relation is a simple 
quadratic equation in w^, and it is straightforward to show that this quantity 
must be real. Stability may then be investigated by passing through the point 



CO 



0. Therefore, marginal instability corresponds to 



and instability is present when 



d[2^ 
d\nR 



= 



dn^ 

dlnR 



< 0. 



(75) 



(76) 



Of course, this condition is nearly universally satisfied for any type of as- 
trophysical disk. This is the simplest manifestation of the magnetorotational 
instability, or MRI. 



4.2 Physical Interpretation 

The equations of motion (j70[) and (|7ip have a very simple mechanical analogue: 
they are exactly the equations of two point masses in orbit around a central 
body, bound together by a spring of frequency kvA (Balbus & Hawley 1992). 

This mechanical analogy immediately offers a physical explanation of the 
MRI (see figure [3]). Two masses are connected by a spring, an inner mass 
and an outer mass mo- The mass orbits faster than the outer mass nio, 
causing the spring to stretch. The spring pulls backwards on nii and forwards 
on nio- The negative torque on vrii causes it to lose angular momentum and 
sink, while the positive torque on mo causes it to gain angular momentum and 
rise. Thus, even though the spring nominally supplies an attractive force, in a 
rotating frame this drives the masses apart. The continued stretching of the 
spring just makes matters worse, and the process runs away. The mass with 
the higher angular momentum obtains yet more, the mass with less angular 
momentum loses what little it has. It is the "same old story" in a dynamical 
context: the rich get richer and the poor get poorer. 

In a protostellar disk, it is the magnetic tension force that plays the role 
of the spring, and the two masses are any two fluid elements tethered by 
the field line. The linear phase of the instability is the exponentially growing 
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Fig. 3. The magnetorotational instability. Magnetic fields in a disk bind fluid el- 
ements precisely as though they were masses in orbit connected by a spring. The 
inner element rrii orbits faster than the outer element rrio, and the spring causes 
a net transfer of angular momentum from rrii to mo. This transfer is unstable, as 
described in the text. The inner mass continues to sink, whereas the outer mass rises 
farther outward. (Figure courtesy of H. Ji.) 

separation of the displaced fluid elements, which is followed by the nonlinear 
mixing of gas parcels from different regions of the disk. The mixing seems to 
lead to something resembling a classical turbulent cascade, though the details 
of this process, with different viscous and resistive dissipation scales, remain 
to be fully understood. 

Notice that angular momentum transport is not something that happens 
as a consequence of the nonlinear development of the MRI, it is the essence 
of the MRI even in its linear phase. The very act of transporting angular 
momentum from the inner to outer fluid elements via a magnetic couple is a 
spontaneously unstable process. 

4.3 General Adiabatic Disturbances 

If /2 is a function only of cylindrical radius R, then for general magnetic field 
geometries, local incompressible WKB disturbances with space-time depen- 
dence 
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5X exp(ifci{-R + ikzZ — wt) 



(77) 



where 5X is the Eulerian perturbation of any flow quantity, satisfy the fol- 
lowing set of linearized dynamical equations: 



iuj5vR — 2fi5vA, + i— 
P 



knSvR + kzSvz = 0, 
B • SB 



6P + 



4n 

k\ k-B 



i— oBr 



— IWOVz + 1 

P 



5P + 



B ■ SB 



An 



4np 
SBr = 0, 
k ■ B 



0, 



6B^ = 0. 



The linearized induction equations are 

- icoSBR - i{k • B)SvR = 0, 



iuSBd, — 5Bt 



i{k • B)5v^ = 0, 



dlni? 

— iuiSBz — i{k • B)6vz = 0. 
Finally, the entropy equation is 

tu'y h Sv;, = 0, 

p oz 



(78) 
(79) 
(80) 
(81) 

(82) 

(83) 
(84) 

(85) 



where 7 is the adiabatic index (not to be confused with the collsion rate). We 

have ignored background radial gradients in the entropy and pressiire, but 
retained their vertical gradients, in accordance with the assumption of a thin 
disk. 

This set of linearized equations leads to the dispersion relation (Balbus & 
Hawley 1991) 



1,2 



kR 



Co^ - An^{k • va) 



0. 



Here 
and 



u>^=uj^-{k-VAf, k-'^kl + kl, 
1 dP d In Pp-'^ 



= — 



dz 



(86) 



(87) 



(88) 



7/9 dz 

N is known as the Brunt- Vaisala frequency. It is the natural frequency at 
which a vertically displaced fluid element would oscillate in the disk due to 

buoyancy forces. In general there is also a contribution due to radial gradi- 
ents as well (Balbus & Hawley 1991), but this usually may be ignored in a 
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rotationally supported (supersonic orbital speed) disk. Without any magnetic 
field, the dispersion relation becomes 

= [k./kfn^ + {kR/kfN"^ (89) 

Since the displacement of the fluid element is incompressible, the wave number 
ratio kz/k is a measure of the radial displacement, whereas fc^/fc is a measure 
of the werticaZ displacement. There is no instability in this case, only a wavelike 
response due to the restoring Corilois and buoyant forces. 

Notice that N vanishes in the disk plane by symmetry, it affects only the 
behavior of disturbances at least one scale height or so above z = 0. Moreover, 
the actual value of N is likely to be determined by radiative diffusion processes 
in the vertical direction. The radiation requires a source, which in disks is the 
turbulence we are trying to explain in the first place! (In protostellar disks, 
external heating of the upper layers is also a source of departure from adiabatic 
behavior.) 

The other limit of interest is k^ = 0, which returns us to the disper- 
sion relation of the previous section. (The expression kvA should of course 
be replaced by k • va-) This is fortunate because it can be shown (Balbus & 
Hawley 1992) that the most rapidly growing modes are displacements in the 
plane of the disk with kji = 0. Since it is the differential rotation that ulti- 
mately destabilizes, it is very sensible that these displacements, which most 
effectively sample the differential rotation, are the most unstable. For a Keple- 
rian rotation profile, it may be directly computed from the dispersion relation 
that the wave number of maximum growth is given by 

k • va = (\/T5/4)i7 (wavenumber of maximum growth), (90) 

or 0.97^2. By comparison, the largest unstable wavelength (i.e., smallest un- 
stable wave number) has a value of 1.73^? for k • va- The maximum growth 
rate is 

|a;| — (3/4)]7 (maximum growth rate) (91) 

This is a very large growth rate, with amplitudes growing at a rate of 
exp(37r/2) ~ 111 per orbit. 



4.4 Angular Momentum Transport 

The azimuthal equation of motion (j45p may be written 



dt 



R pVV^ - pVAVA<t, 



P+ — 
Svr 



0, 



(92) 



where va4> is the azimuthal component of the Alfven velocity. This is an equa- 
tion for angular momentum conservation, with angular momentum density 
pRvcj), and an azimuthally averaged flux of 
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Fj = {R {pvv^ - pvAVA<t>)) (93) 

In protostellar disks, we are interested in the radial transport of angular mo- 
mentum, 

FjB. = Rp{vRV^ ~ varva4>) (94) 

It is most instructive to begin with the simple case of linear instability in 
a uniform, vertical magnetic field. What is the lowest order flux that results 
from these exponentially growing modes? In equilibrium, vr and the Alfven 
velocities vanish. The linear pertubation of the radial velocity at a given spa- 
tial location (a so-called Eulerian perturbation) will be denoted 5vr. Let the 
radial displacement of a fluid element from its equilibrium location be given 

by 

= ae'^* cos(fcz) (95) 

where 7 and k correspond to the maximum growth rate and its associated 
wave number, and a is a slowly varying amplitude. Then 

SvR. = Ci? = 1S.R (96) 

The azimuthal velocity consists of an unperturbed Keplerian velocity 
plus a linear perturbation 5v^. The product of 5vr and J7 contributes zero 
when a height integration is performed because of the cosine factor, so we must 
consider the direct product 5vr x 5v^ to obtain the lowest order contribution. 
In calculating ^z;^, note that 

The subtraction is needed to eliminate the change in velocity a displaced fluid 
element would make even if there were no physical change in the rotation 
velocity at the new radial location: the actual change in velocity 5v^ is due to 
the change in the displaced fluid element makes in excess of ^Rdf2/dhi R. 
Equation ([70|l then gives 

^^^ = ^^-^-dh^=(^-^j^« ^'^^ 
and therefore the velocity contribution to the angular momentum flux is 

{5vR5v^) = ^{l^ + ev\){eR) (99) 

Not surprisingly, it is positive (outward). 

To calculate the magnetic field correlation, we use equation ((S7)) to go 
between magnetic field fluctuations and displacements. Then, using (j7ip . we 
find 

U = - 2 2 (100) 
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from which it foUows, among other things, that and dvc/, have opposite 
signs. Hence 



(CI) 



Therefore 



7 

2Q 



7^ + fc2x;2 



7^ + k^v\ 



(101) 



(102) 



Using the dispersion relation (j74p to simphfy this a bit, our expression for the 
radial angular momentum flux is 



(103) 



7^ + k'^v\ 

Finally, integrating over height and defining the effective surface density S by 



(104) 



the angular momentum flux for the most rapidly growing mode is found to be 



FjR 



(105) 



This result bears some commentary. First, the radial flux is always positive 
for any unstable mode, not just the most rapidly growing. This, as we shall 
see, is a direct consequence of energy conservation: energy is extracted from 
the differential rotation to the fluctuations only if the radial angular momen- 
tum flux is positive. Next, note that there is an outward angular momentum 
flux only to the extent that there is a correlation in the velocity fluctuations. 
This is also true in turbulent flows, and everything that we have done in this 
section goes through in much the same way if the fluctuations arc turbulent as 
opposed to wavelike. The important physical point is that an instability that 
leads to turbulence need not lead to enhanced angular momentum transport. 
Only turbulence with strongly correlated velocity fields does this. Strong cor- 
relations, in turn, are necessary to extract energy from differential rotation. 
This is indeed the free energy source of shear-driven turbulence, but not any 
form of turbulence. 

This point was made in Stone & Balbus (1996), in which the angular 
momentum transport resulting from convective instability was studied. The 
angular momentum was tiny in magnitude and inwardly directed: it had the 
"wrong" sign. Disk intabilities based on adverse thermal gradients do not, 
in general, lead to systematically large outward transfers of angular momen- 
tum. When it comes to angular momentum transport, not all turbulence is 
equivalent. 



24 Steven A. Balbus 



4.5 Diffusion of a Scalar 

A classical problem in protostellar disk theory is to understand how dust 
particles are mixed with the gas. The fact that the MRI leads to vigorous 
radial angular momentum transport suggests that other quantities may be 
transported as well. In this section, we will give an argument that shows 
that the radial diffusion of angular momentum is indeed closely related to the 
radial transport of a conserved scalar quantity, say Q. Assume that Q is a 
fluid element label and satisfies an equation of the form 

d \ 

+ t;.V = (106) 



Then, 



dt 

Using mass conservation, this implies 



dt 



^ + V-{vQ)-QV-v^Q. (107) 



dt ^ ^' ^ Dt 

where 



' v.(.g) + g^^o, (108) 



D_ fd_ 

[di. 

Simplifying the equation for Q, we obtain 

d{pQ) 



v-V (109) 



, V-{pvQ) = 0, (110) 

which looks very much like the angular momentum conservation equation (|92p . 
In a similar way, there will be a turbulent Q-flux given by 

Fq = p{SvSQ) (111) 

Since Q is conserved as we follow a fluid element, the so-called Lagrangian 
perturbation AQ, defined by 

AQ = SQ + ^-VQ (112) 

must vanish. This is because AQ is constructed to be the change in Q following 
a fluid element as it is displaced by a small distance ^, and such changes must 
vanish since by definition Q is does not change along fluid element paths. 
Hence, our expression for the ith component of Fq is 

Fq, = -p{Sv,^^)djQ (113) 

This defines the diffusion tensor 1?^ , 



(114) 
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In a height-integrated calculation, it is the RR component of V that is of 
importance. In the linear regime 6vfj — j^r, so that 

Vrr = 7(efl) (115) 

In this case, equation (|103p gives a relationship between the flux of angular 
momentum and the diffusion coefficient of a passive scalar. 



7^ + fc^w^ 



(116) 



For the fastest growing linear mode, this gives Fjh = 2fiT>iiB_. 

There is a sense in which the linear theory might find its way into nonlinear 
turbulent diffusion. In simulations the MRI appears to locally stretch field 
lines and exponentiate velocity growth over a limited duration of time, before 
one fluid element becomes mixed with another and the process starts anew. 
If such a picture is reasonably accurate, then Fjr ^ i^I?_R_R may well be valid 
in turbulent flow. (It is also, of course, what we might expect on the basis of 
simple dimensional analysis alone!) The important astrophysical point is that 
if protostellar disks are MHD turbulent, they ought to be well-mixed. 



5 Energetics of MHD Turbulence 
5.1 Hydrodynamic Considerations 

The quantity 

{pVV^, - pVAVA^) (117) 

plays two conceptually different roles in the theory of accretion disks. We have 
seen that it is intimately linked to the direct transport of angular momentum, 
and in this section, we shall study in detail how it extracts free energy from 
the large scale differential rotation. 

Let us begin with the relatively simple case of an adiabatic nonmagnetized 
gas. The azimuthal velocity is decomposed into a time-steady large rotational 
component Rf2 plus a fluctuating component u^: 

v^ = RQ + u(j) (118) 

we will assume neither than u(j) is small compared with RQ nor that the mean 
value of Uff, vanishes, though in practice both might well be the case. We will 
write ur and Uz for the radial and vertical velocity components respectively, 
and V for the full velocity vector. Once again, we think of these as fluctuations, 
but our treatment will in fact be exact. 

The adiabatic radial equation of motion is 
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Multiplying by Ufj and regrouping: 
'd 



dt 



where 



R 



dP 
dR 



eff 



OR 



eff 



Exactly the same manipulations with the 



d_ 

dt 



— +vV 



PU^UR 



sQ'^{s) ds 

and z equations produce 
dP 



2n 



} 

R 
dP 



Rd(t> 



P- 



eff 



dz 



(120) 



(121) 



(122) 



(123) 



where we have assumed that <P is independent of 0. Adding the three dynam- 
ical equations and using mass conservation leads to 



d_ 

dt 



1 



-pu 



V-{pu^v/2) + puRU^ 



df} 
d\nR 



u ■ VP - pu ■ V^, 



eff 



(124) 



where = 



-u^. Yet another use of mass conservation and a regouping 



of the pressure term gives us 



d /I 



2^" 



P'^eff 



V- ( V{^pu'^ + p^eff) + Puj = PV-U- pURU^-^^^ 



d\nR 
(125) 

We have already at hand the main structure of our energy equation. The left 
side is in conservation form with a well defined energy density and energy 
flux, with the fluctuations isolated. (Notice that the azimuthal average of the 
energy equation, which we ultimately will be working with, has no rotational 
terms in the energy flux divergence.) Sources of energy fluctuations are work 
done by pressure (which may cause heating or cooling), and the all important 
R(p stress coupling to the differential rotation. 



5.2 The EfTects of Magnetic Fields 

When magnetic fields are included, everything goes through as before, except 
that the right side of equation (|124p contains the terms 

- u ■ VPtot - pu ■ V^ef f + "'[^ • (126) 

47r 

The final magnetic term in the above may be written in an index notation as 

u^BABj^ (127) 
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where k take the values x, y, z, and repeated indices are summed over, di 
denotes the the partial derivative with respect to the ith spatial variable. 
To make further progress, we need the induction equation: 

^^+vV^B = -BV-v + {B-V)v. (128) 

Take the dot product of this with B and write the last term in component 
form: ^ ^ 

|_ + vv) ^ = ~B^V-v + BjBAvj (129) 



dt 



Now, 



BjBAvj = d.i[B^BjVj)-v.jd^{BjBi) = V-{v-BB)-R[2[B-V B]^-UjB,d^Bj 

(130) 

where we have used diBi = 0. Note that with the last term, we make contact 
with the energy equation terms p26p . The subscript (j) on the penultimate 
term denotes a vector component: 

RQ[B-VB]^ = n[B-V{RB^)] = V-{RnBB^) - B^Br-^ (131) 

a In rt 

Therefore, the right side of equation p29p becomes 

- B^V-v + V-[{u-B)B] + B^Br-^ - UjBAB, (132) 

a m rt 

Combining this result with the left side of the equation p29p then, after some 
cancellation, given us an expression for ujBidiBj: of terms, 

-u.BAB, - (I + ^ + V. (u^ - in . B)B 

+ fv..-5,i?«^ (133) 

Armed with this result, we return to (|126p and make a substitution for the 
final term. The resulting energy equation can be simplified to 

9 2 ^ V, r /^l 2 ^ \ n Bx(uxB) 



V -pu' + p-P.ff + Pu 



2'^ ' / 47r 



+ V- {e^R[2— \ = PV-u - p{uiiu^- VARVA^)-^^^^ (134) 

The pressure term on the right side of the above equation can be eliminated 
using the thermal energy equation 

^^+v-v\- = -PV-u- pC (135) 
.dt J p 
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(since V-t) = V-tt). We have introduced the radiative energy loss term per 
unit volume pC Carrying through the elimination of the pressure term and 
averaging over (j) lead to the total energy equation: 



^ + V-J^ = -p{uRU^-VARVA4,)-^l^^- pC (0 averaged) (136) 
where £ is the energy density in fluctuations 



1 B'^ 3P 

£ ^^pu^ + P'Peff + -^ + ^, (137) 

and ^ is the corresponding flux, 

_ /I o \ 5P B X (u X B) , , 

J^ = ui^-pu' + p<P,ffj+—u+ (138) 

We have completely ignored dissipation effects (viscosity and resistivity). 
What effect would this have on our final equation (|136p ? The answer is es- 
sentially none. While it is true that new dissipation terms would appear on 
the right side of equation (|134p . they would also appear with the opposite 
sign in equation (|135p . They would completely cancel on the right side of the 
final equation (|136p : dissipation is not a loss of energy, but a conversion of 
mechanical or magnetic field energy into heat. Additional small energy flux 
terms would appear (within the divergence operator) in connection viscosity 
and resistivity, but these are generally negligible compared with the dynam- 
ical terms. The essential physical point is that total energy is conserved in 
the presence of dissipation even if mechanical energy is lost. Only radiative 
processes can remove energy from the disk. 

We interpret equation p36p as saying that energy is exchanged with the 
differential rotation at a volumetric rate —Tj^^dQ/dlnR^ where 

Tr^ = p{uRU^ - VARVA4>) (139) 

is the dominant component of the stress tensor. The energy made available 
from the differential rotation may in principle remain in velocity fluctuations 
in the form of waves, but in a thin, cool disk it is more likely that this energy 
will be locally dissipated, and subsequently radiated. This is because if the 
energy flux varies radially over a scale of order R itself, then the T^^ source 
term on the right will be larger than any term on the left (by an amount of 
order u/{RQ)). The stress can only be counterbalanced by the radiative loss 
term: 

- TR^dn/dlnR = pC (140) 

Therefore, although —Tii^df2/dhiR is itself a nondissipative source term, in 
these so-called local models, it works out to be the rate at which energy is 
dissipated (and ultimately lost by radiation) as well. Bear in mind that the 
identification of the stress term with dissipation follows from a reasonable 
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assumption about the behavior of thin disks: local dissipation of the fluctu- 
ations. Dissipation is in no way a fundamental property of the left side of 
equation (|140p . which under different conditions could just as well be an en- 
ergy source or reserve for purely adiabatic waves. There are many papers in 
which this point is misunderstood: caveat lector. 



6 Resistive and Hall terms 
6.1 Local Dispersion Relation 

Extensive regions in real protostellar disks are in a regime far from ideal MHD, 
in which Ohmic decay and Hall electromotive forces are important (e.g. figure 
1). It is important to understand how the MRI behaves under these conditions. 

The fastest growing modes, as usual, correspond to axial wavenumbers 
k = ksz ■ Fortunately, these are also the simplest to treat analytically, and we 
will limit our discussion to this class of disturbance. 

The linearly perturbed Hall term in the induction equation is to leading 
WKB order 



- Vx 



dJ X B 


= Vx 
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since 6 J = (c/47r)(ifc X 6B). Taking the curl and simplifying allows us to 
write the right side as 

^ {k-B){kxSB). (142) 



47ren,_ 

(In a fully ionized disk, the presence of Uf, in the denominator renders these 
terms negligbily small.) For our particular problem, k is axial and SB lies in 
the disk plane, and the above reduces to 

^—{SBRe^-SB^en). (143) 

Including the effects of resistivity, the linearized induction equations (|82ll83p 
become: 

ck'^B 

{-iuj + rik^)SBB + SBs - i{k ■ B)Svb = 0, (144) 

(_,^ + rjk')SB^ - SBn — — + ) - z{k ■ B)Svr = 0, (145) 

and the z equation is not needed (trivially satisfied). The linearized dynamical 
and energy equations from section (4.3) remain unchanged. 

The Hall electromotive force terms introduce a phase shift in the induction 
that leads to a circularly polarized component of what would otherwise be 
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an Alfven or slow mode. When an instability is present due to differential 
rotation, this phase shift can either destabilize or stabihze depending upon 
the sign of B^. 

To see this more quantitatively, we must study the dispersion relation 

that follows from our linearized set of equations. To this end, we define a Hall 
parameter vh which has the dimensions of a velocity (Balbus & Terquem 
2001), 

v% = TT^- 146 

Notice that if we take the sense of rotation to be positive, that v\ can have 
either sign depending upon whether is positive or negative. With u = —iu, 
our dispersion relation is 

+ 2r]k'^ + 620^ + 2r]k'^{K^ + k'^v\)a + Co = 0, (147) 

where the constants C2 and Cq are given by 

C. = .^ + 2.M + .^^^ + #(^+^X) (148) 
Co = {v\ + ^) + + k'v%^ + ^Wk' (149) 



6.2 Stability 

A sufficient condition for instability is that Co < 0. When compared with the 
purely Ohmic requirement 

k\\ ( + + i^Wk^ < (150) 

\ dhiK J 

The Hall condition Co < is much more easily satisfied, since can change 
sign. Perhaps most striking is the possibility that the intuitive result that 
large wavenumbers are stabilized can be lost when v'jj < 0, and that all 
wavenumbers can be unstable, even in the presence of finite resistivity. (See 
figure [4].) This may be shown straightforwardly from the definition of Co; the 
reader may wish to consult Balbus & Terquem (2001) for explicit details. In 
addition to requiring counteraligned magnetic and rotation axes, this state 
also demands a magnitude ratio of v\ to v'jj lies between about 2 and 4 (for a 
Kcplcrian disk). When the field is aligned with the rotation axis, the range of 
unstable wavenumbers is more restricted than that found for the ideal-MHD 
MRI. 
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Fig. 4. Range of unstable wavenumbers for a representative disk in the Hall-Olim 
regime. Tlie black diamond corresponds to the MRI without the Hall effect. The 
resistivity is chosen to be K^rf /v\ = 0.35. Prom Balbus & Terquem (2001). 



6.3 Numerical Simulations of the Hall-Ohm-MRI. 

What does the preceding section imply for the MHD transport properties of 

protostellar disks? The first local simulatons of the nonlinear development of 
the MRI in the presence of Ohmic resistivity are those of Fleming, Stone, & 
Hawley (2000). These authors found the important result that the nonlinear 

turbulent state is easier to maintain when a mean vertical magnetic field is 
present. This is sometimes misunderstood as a requirement that a disk must 
have a global vertical magnetic field to maintain MRI turbulence, but in fact 
all that the very local calculations of Fleming ct al. (2000) really tell us is 
that on this scale, working with a precisely zero mean field is probably too 
restrictive an approximation for a real disk. This study found a difference of 
two orders of magnitude in the critical magnetic Reynolds number (defined 
here as the product of the isothermal times the box size divided by the resis- 
tivity) needed to sustain turbulence, depending upon whether a vertical field 
was present or not, the vertical field runs being the easier to sustain. It is now 
understood that very high resolution grids are needed for the local study of 
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the MRI in a shearing box with zero mean field, so these early results should 
not be used quantitatively. With this precaution understood, it is interesting 
to note that for the vertical field runs, the Fleming et al. finding seems to be 
not very different from the L = 1 criterion discussed earlier. 

This nature of the nonlinear behavior of a Hall fluid vi^as examined in a 
numerical investigation by Sano & Stone (2002a, b). These authors carried 
out an extensive study of the local properties of Hall-MRI turbulence in the 
shearing box formalism. The results arc somewhat surprising. 

In linear theory, counteraligned magnetic and rotation axes produce a very 
broad wavenumber spectrum of instability, and should be more unstable than 
the aligned case. What Sano and Stone found, however, was just the oppo- 
site: the aligned case showed greater levels of field coherence and higher rms 
fluctuations than the counteraligned case. The reason for this is interesting, 
and illustrative of the hazards of extrapolating directly from linear theory. 

In the case of counteraligned axes, the extended wavenumber of spectrum 
of linear MRI instability results, in its nonlinear resolution, in a highly ef- 
ficient turbulent cascade to smaller and smaller scales. In a finite difference 
numerical simulation, this cascade terminates at the grid scale where energy is 
ultimately lost. By contrast, the aligned case shows instability only at longer 
wavelengths (small wavenumbers) . At larger wavenumbers, the response of 
the fiuid is wavelike, and this is likely to make a local small scale turbulent 
cascade considerably less efficient. The result: a less lossy system and greater 
large scale field coherence. 

Finally, when simulations were done with half the field aligned with the 
rotation axis and half counteraligned, the Hall effect tended to wash out, and 
the results were similar to resistive MHD turbulence without the additional 
electromotive forces. The criterion for the onset of turbulence depended only 
the relative strength of the resistivity, as measured by the dimensionless num- 
ber v\/rifi. This so-called Lundquist nvmber must bo greater than imity or 
MRI-induced turbulence is not maintained, whether Hall electromotive forces 
are present or not, according to Sano & Stone. The presence of the Hall terms 
in the induction equation do, however, promote a slightly elevated rms level 
when turbulence can be maintained. 

The Hall numerical studies are fascinating and suggestive, but caution is 
necessary before taking the final leap from simulation to real disks. Note, for 
example, the role played by the numerical grid scale losses in our discussion 
of the nonlinear resolution of the counteraligned case. Nature works without 
a grid, however, and it is not clear what scale ultimately intervenes at high 
wavenumbers in Hall protostellar disks. Moreover, it is now understood that 
there are important resolution questions for local shearing box simulations 
with no mean field, that have only recently been brought to light (Fromang 
et al. 2007; Lesur & Longaretti 2007). The important pioneering studies of 
Sano & Stone (2002a, b) merit a thoroughgoing follow-up on the much larger 
numerical grids that are currently available to numericists. 
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7 Alpha models of protostellar disks 

7.1 Basic principles 

We are now in a position to assemble a basic "a model" of a protostellar disk. 
An important assumption in such models is that angular momentum and me- 
chanical energy are transported radially through the disk. We therefore work 
with height and azimuthally integrated quantities. Under these conditions, 
Balbus & Papaloizou (1999) have shown how the fundamental equations of 
MHD lead directly to the a formalism. 

The velocities consist of a mean part plus a fluctuation of zero mean. For 
the azimuthal velocity the mean will be taken to be the Keplerian velocity 

W = ^ (151) 

where M is the central mass and G is the Newtonian constant. For the radial 
velocity fij, the mean motion corresponds to the very slow inward accretion 
drift, vi. The azimuthal fluctuation velocity is much less than _Ri7, whereas 
the radial fluctuation velocity ur is much larger than the inward drift. These 
claims will shortly be quantified. The magnetic field Alfven velocities are 
assumed to be of the same order as the u velocities, namely less than or of 
order the isothermal thermal sound speed cg. To summarize, 

v-i. <C «, VA, cs «: RO. (152) 

Under steady-state conditions, the local disk accretion rate has a time 
averaged value of 

= 2nR{pv2 + 5p5uR) (153) 
The height integrated form of this is — M, the mass accretion rate: 

IJ,dz = -M, (154) 



I' 



defined so that M is a positive constant. A useful average of V2 is the drift 
velocity defined by 

where S is the disk surface density. 

In steady-state, the height-integrated average radial angular momentum 
flux is proportional to 1/R. In other words, 

- MR'^Q + 2'kSRWr4, = C/R (156) 

where Wr^ is deflned by 

{pURU^ - pVARVA<t>) dz = SWr^, (157) 
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and C is a constant to be determined. Traditionally, this has been done by 
asserting that W^^ is proportional to the radial gradient of SI, and at some 
point before the surface of the star is reached this must vanish. C is then 
determined at this point (Pringle 1981). 

This approach now seems dated, and especially inappropriate for a proto- 
stellar disk in which the magnetic interactions between the disk and star may 
become quite complex at small radii. Instead, let us note that at small R, C / R 
must be very small if SWr^ does not blow up, which seems quite reasonable. 
If C/i? is very small at the inner edge of the disk, it is clearly negligible in 
the body of the disk, and we shall set C — with the understanding that our 
solution should not be taken to the stellar surface. Then, (|156p reduces to 

Vd = -WH4/Rn (158) 

showing that Vd (and V2 which is of the same order) is very small compared 
with the u velocities. Alternatively, 

SWr4, = MRf2/27r. (159) 

This is a particularly useful result because Wu^ appears in the energy balance 
equation (|140p : 

- ^^-^d^ - j^-^'^ ^jpCd^^ 2.Tt,, (160) 

In the last equality, we have equated the radiated energy per unit surface 
of the disk to twice the blackbody emissivity, the factor of two represents 
two radiating surfaces. Here cr is the Stefan-Boltzmann constant and T^^fj is 
the effective blackbody surface temperature. Combining (|159p and (|160p then 
gives us 

4 3GMM 

relating the potentially observable disk surface temperature to the central 
mass, accretion rate, and radial location. The unknown turbulent stress pa- 
rameter Wi?^ has conveniently vanished! 

Equation (|16ip is in a restricted sense "exact." It is based on the assump- 
tion of time-steady conditions and thermal radiation, but it is independent of 
the explicit nature of the turbulence, as long as it is local. If these conditions 
are all met, equation (|16ip is a simple matter of energy conservation. To go 
beyond this result, something has to be said directly about W/j^. In a disk 
theory, that condition is 

Wr^ = acl (162) 

where a is assumed to be constant, but otherwise unknown. In contrast to the 
simple and plausible assumption that turbulence is locally dissipated, this "a 
assumption" is far from obvious, and in a strict mathematical sense almost 
certainly wrong. 
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The original justification for this form of the stress tensor was based on 
the notion of hydrodynamical turbulence and the idea that the velocity fluc- 
tuations would be restricted to some fraction of the sound speed (fluctuations 
in excess would cause shocks). Shakura & Sunyaev (1973) consider magnetic 
stresses, however, and argue that they fit within the a formalism as well since 
Alfven velocities in excess of cs are also dynamically unlikely. 

The real problem with the prescription (|162p is that turbulence is just too 
complicated. Not only are long term averages very hard to define (a problem 
even for a result like equation |161j ). it is also entirely possible that a could 
vary in a complex nonsystematic way by an order of magnitude or more from 
one part of the disk to another. The primary justification for (|161[) is that 
many scaling results are very insensitive to a. On dimensional grounds cs is 
is certainly an important local characteristic velocity, but it is not yet clear 
under what conditions a "background" magnetic field might also be providing 
a mean Alfven velocity that limits or guides the turbulence. 

Continuing with our a model, to link the surface temperature T^ff with 
the midplane temperature T (used in quantities like cs) we need to intro- 
duce explicit vertical structure into the problem. The hydrostatic equilibrium 
equation is 

dP GMpz 2 

This may be solved in conjunction with a detailed energy equation, but it 
seems best for illustrative purposes to note that this equation serves to provide 
a decent and very simple estimate of the disk half thickness, H — cs/fi- 
The energy equation for simple vertical radiative diffusion defines the local 
radiative energy flux F-y as: 



Fi = — — - (164) 



4(7 dT-* 

The "optical depth" r is given by 

dr = —pudz (165) 

where n is the opacity of the disk (units: cross sectional area per unit mass). 
In the simplest possible model, is a constant given by uT^^^, and the 
temperature at the midplane T is then 

= ^ (166) 

4 

where r is the optical depth from the outer disk surface to the midplane. (In 
this simplest of all possible models, we set r — kH, where k is evaluated at 
the midplane temperature.) 

The surface temperature of a disk with a mass accretion rate of 10~^Mq 
per year around a lA'/© star is SbR'^'^'^^K, where Rau is the radius in astro- 
nomical unit. The midplane temperature T is a factor ~ t^^^, larger, typically 
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a factor of 5 or so larger. At what value of T would we expect the ionization 
fraction to reach the critical value of 10~^^ we found earlier, corresponding to 
a Lundquist number of unity? Put in somewhat different terms, is our model 
of self-sustaining MHD turbulence self-consistent? 

To answer this, we need to address the physics of thermal ionization. In 
the low ionization regime we are working, thermal electrons are supplied by 
trace alkali elements, notably potassium, with an ionization potential of only 
4.341 eV. Even this modest value corresponds to an effective temperature of 
50,375 K, well above the range of 10^ to lO'^K we expects, and potassium will 
barely be ionized. 

The equation governing the ionization fraction x is the Saha equatiorQ. In 
this barely ionized regime, it may be written (Stone et al. 2000) 

^asi j T^'^ exp(-50375/r) (167) 

where a is the abundance of potassium ( 10^^), n is the ratio of the dominant 
H2 molecules, and s is a ratio of statistical weights, expected to be near unity. 
We may rewrite this equation as 

--^ 

where 

X = x^T-^l'^injlA X 10^^) (169) 

With x — 10^^^ and the final density factor anything between 0.01 and 1 
(O.IM0 of gas spread out over a region of 10 AU yields a density of about 
5 X 10^'^ per cc) gives values for T close to 10"^K, which we will adopt as a 
working number. 

Therefore, a surface temperature of about 200K is required to attain a 
midplane temperature of 10'^, and thereby ensure a critical level of ionization. 
This is the key issue for MHD theories of protostellar disks: beyond 0.2—0.3^1 U, 
or ~ 3 X 10^^ cm, the heat generated by the dissipation associated with local 
turbulence is insufficient to keep the the disk magnetically well- coupled. Where 
and how do protostellar disks maintain good magnetic coupling? 



8 Ionization Models of Protostellar Disks 
8.1 Layered accretion 

The above considerations indicate that dissipative heating from MHD turbu- 
lence self-consistently generates enough heat to maintain the minimum ther- 

^ The discussion presented here presumes that the ionization reaction rates are 
sufficiently rapid to keep up with the dynamical time scales of any turbulence 
present. At low ionization fractions this breaks down, and a single temperature 
may not be enough to characterize the ionization state (Pneuman & Mitchell 
1965). 
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mal levels of requisite ionization only with a few 0.1 AU. Nonthermal sources 

of ionization arc therefore of great interest to protostellar disk theorists. The 
principle ionizing agents that have been studied are cosmic rays, X-rays, and 
radioactivity. 

Gammie (1996) argued that the low density extended vertical layers of a 
protostellar disk would be exposed to an ionizing flux of interstellar cosmic 
rays. Just as in models of molecular cloud ionization, cosmic rays would main- 
tain a minimal level of ionization in the upper disk layers. The range of the 
low energy galactic cosmic rays is about 100 g/cm^. This is much less than 
the disk column density at ~ 1 AU in generic solar nebula models, but can 
easily exceed the disk column at larger distances. If the level of ionization is 
high enough-and the Alfven velocity va is small enough — the gas within the 
range of the cosmic rays will be MRI active. (The Alfven velocity cannot be 
too large if disturbances are not to be stablized by magnetic tension.) If these 
criteria are met, Gammie argued that turbulent accretion would proceed in 
the outer layers of protostellar disks, but that the midplane regions would 
remain laminar in effect, a "dead zone." This is the basis of the concept of 
layered accretion, which has become an important idea in protostellar disk 
modeling (Sano et al. 2000). 

Gammie's original construction was based on the assumptions that the 
accretion would occur in a layer of fixed column, and that a is constant. 
Taken together, these assumptions preclude a steady solution; the mass flux 
rate is not independent of position. Instead, matter is deposited from the outer 
regions into the disk's inner regions where it builds up. At some point in this 
scenario the disk become gravitationally unstable, and it was speculated that 
an accretion outburst might occur, which was tenatively identified with FU 
Orionis behavior. On the other hand, there is no compelling argument (beyond 
mathematical convenience) to prevent the a parameter from adjusting with 
position, if this allows a relaxation to a time steady solution. The overarching 
layered accretion picture would, however, remain intact with this modification. 

YSO's are almost universally X-ray sources. Glassgold et al. (2000) noted 
that X-rays are potentially a far more powerful ionization source, even if atten- 
uated, than galactic cosmic rays. These authors were thus the first to draw the 
link between X-ray observations and MRI activity of accretion disks. In the 
Glassgold et al. formulation, X-rays from a locally extended corona centered 
on the YSO irradiate the disk, and depending upon the chosen parameters, 
could provide the requisite ionization levels to lessen the extent of the dead 
zone or eliminate it altogether. 

Whether the dead zone persists in the presence of X-ray irradiation or not 
depends, among other things, upon the model adopted for the disk structure. 
Generally, the so-called minimum mass solar nebula model (Weidenschilling 
1977) is used. This is a reconstruction of the surface density distribution of 
the sun's protostellar disk based on current planetary masses and compo- 
sitions, and results in an R^^/^ radial distribution. Fromang, Terquem, & 
Balbus (2002) suggested that ionization fractions should be calculated self- 
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consistently within the framework of accretion disk theory, which in general 
leads to a more shallow dependence of smfacc density with radius. By intro- 
ducing accretion parameters whose values are free (a and the mass accretion 
rate M), the range of parameter space increases, and the presence and extent 
of a dead zone becomes yet more model dependent. 

The ionization fraction of a protostellar disk is a chemical problem, and in 
principle can involve a very complex, uncertain, molecular reaction network. 
Ilgner & Nelson (2004a) investigated the consequences for the dead zone of 
a considerably richer chemical network. While the quantitative details were 
sensitive to the chemistry, the qualitative structure was not. A dead zone is 
likely to be present in any plausible chemical network that has been studied up 
to the present, but there are conditions in which its extent can be very small or 
possibly even zero. Ilgner & Nelson (2004b), for example, made the interesting 
comment that flaring activity by the central star can have a significant effect 
on the extent of the dead zone. 

8.2 Activity in the Dead Zone 

Just because the dead zone is unable to host the MRI does not mean that 
it is well and truly dead. Fleming and Stone (2003) carried out numerical 
simulations in which the magnetically active upper layers of the disk cou- 
pled dynamically to the magnetically inert dead zone, resulting in a small but 
significant Reynolds stress, even though the Maxwell stress was zero. More re- 
cently. Turner & Sano (2008) suggested that the high resistivity characteristic 
of the dead zone is also a means for magnetic field to diffuse into this region 
from the active layers. Moreover, these authors point out, an active MRI re- 
gion is not an absolute prerequisite for accretion: direct magnetic torques on 
much larger scales would also serve, and since they do not dissipate energy in 
a turbulent cascade, are much less costly to maintain. 

Terquem (2008) constructed global models of protostellar disks based on 
the a prescription. The presence of a dead zone was surprisingly nondisruptive, 
provided that a was not too small. Steady solutions were found for values of a 
in the dead zone as small as 10~^ times the active zone value. In these models, 
the dead zone was thicker and more massive than its surroundings, but because 
of the relatively weak a scalings, by less than an order of magnitude. These 
results, taken as a whole, suggest that an embedded dead zone in a protostellar 
disk is not necessarily disruptive to the accretion process. 

8.3 Dust 

In all of our discussions, we have been most negligent by not mentioning the 
effects of small dust grains (e.g. Sano et al. 2000). Let us see crudely why this 
is so. 

Consider a mass M of protostellar disk gas, of which lO^^Af is in the form 
of spherical dust grains of radius r^. If the dust grains have a density of 3 gm 
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cm~^, there are a total of 

= (10^V47r) (M/r3) 

grains in a volume V. The grains present a geometric cross section of ad = 
7rr^ (we ignore here the enhancement due to the induced charge [Drain & 

Sutin 1987]), and the electrons have an average radial velocity of some = 
1.6(fcr/me)^^^. The total dust recombination rate per unit volume is 

{Na/V)n,aav, (170) 

where rig is the number density of electrons. This should be compared with a 
typical dielectronic recombination rate of P = 8.7 x 10~^ cm^ s~^(Glassgold, 
Lucas, & Omont 1986; Gammie 1996). If x is ionization fraction, then the 
ratio of dielectronic recombination to dust recombination is 

NdadVe ^ V 10-13 J T 

where is in cm and T in K. For small dust grains (~ 10~^cm) this ratio 
is <C 1 and dust recombination is overwhelmingly important. But the relative 
importance of the grains diminishes as they grow in size, especially in the 
cooler portions of the disk. Sano et al. (2000) concluded that in their model 
(based on cosmic ray ionization) a typical protostellar disk would be MRI 
stable inside of about 20 AU if small grain are present, except for the innermost 
regions which would be thermally ionized. 

The detailed effects of dust grains have been considered by many authors; a 
very good review and list of references is given in the recent paper of Salmeron 
& Wardle (2008). These authors have considered the vertical structure of the 
magnetic coupling in a minimum mass solar nebula model, including the ef- 
fects of Hall electromotive fores and ambipolar diffusion. They employed a 
sophisticated chemical network (Nishi et al. 1991), and small dust grains 
ranging in radius from 0.1 to 3 microns. At the fiducial locations of 5 and 
10 AU, Salmeron & Wardle found good magnetic coupling over an impressive 
range of magnetic field strengths provided for the 3 micron grains, a consider- 
able imporovement from our naive estimate above. For the 0.1 micron grains, 
the magnetic coupling dropped sharply (as did MRI growth rates), and was 
restricted to higher elevations above the disk plane. 

9 Summary 

Protostellar disks are gaseous systems dynamically dominated by Keplerian 
rotation. These disks are accreting onto the central protostar, so there must 
be a source of enhanced angular momentum transport present. In the early 
stages of the disk's life, this enhanced transport may well be due to self- 
gravity, with density waves largely responsible for moving angular momentum 
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outwards. Once the disk becomes observable, its mass is below the minimum 
needed to sustain self- gravitational spiral wave transport. The only established 
mechanism able to sustain high levels of angular momentum transport is MHD 
turbulence produced by the magnetorotational instability, or MRI. 

The dynamical effects of magnetic fields on protostcllar disks depends cru- 
cially on the degree of ionization, more specifically on the electron fraction, 
that is present. This means that there is a direct link between the gross dy- 
namical behavior of a protostcllar disk and its detailed chemical profile. In 
principle, the full multifluid nature of the disk gas — neutrals, ions, elec- 
trons, and dust grains — must be grappled with at some level to elicit and 
understand the disk strucutrc. Fortunately, only a very small ionization level 
is needed to couple the charged and neutral components, leaving (in essence), 
a single magnetized fluid. An electron fraction of 10~^^ is a typical fiducial 
number for the threshold of magnetic coupling at 1 AU in the T Tauri phase 
of the solar nebula: roughly one electron per cubic millimeter of disk gas! 
Unfortunately, however, on scales of lO's of AU, near the dense midplane 
the ionization of protostcllar disks may not even rise to this minimal level of 
ionization. At the boundary between good and poor magnetic coupling, the 
MHD processes are complex, not well-understood, and the domain of ongoing 
inquiry. 

In regions of the disk where the magnetic coupling is sufficient, the com- 
bination of a weak (subthermal) magnetic field and Keplerian rotation leads 
to the MRI. In ideal MHD, the magnetic field behaves as though it were 
frozen into the conducting gas, and the presence of differential rotation pro- 
duces azimuthal magnetic field from a radial magnetic field. But this is not 
all that happens. If one tries to simulate this simple process on a computer, 
the laminar fiow breaks down into a turbulent mess, even if the field is very 
weak. 

The reason for this is due to another classical property of magnetic fields, 
that the force exerted by the lines of force on the background gaseous fluid 
is the sum of a pressure-like term and a tension-like term. It is the latter 
tension-like term that is important for an understanding of the MRI. When 
two nearby fluid elements are moved apart, even if only because of a random 
perturbation, the magnetic tension force acts precisely like a spring coupling 
the two masses. The fluid clement on the inside rotates faster than the element 
on the outside, and tries to speed it up. The outer element, in acquiring 
angular momentum from the inner element, finds itself too well-endowed, and 
spirals outward toward a higher orbit where its excess angular momentum 
can be accommodated. On the other hand, the inner element, having lost 
angular momentum, finds itself at a deficit and must drop to a lower angular 
momentum orbit. This separation stretches the field lines that couple the 
elements, the magnetic tension goes up, and the process runs away. This is 
the MRI. Notice that angular momentum transport is at the very core of 
the linear instability, rather than the result of some sort of nonlinear mixing 
process. 
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But the MRI does, in fact, lead to rapid turbulent mixing of the disk 

gas, as outwardly moving and inwardly moving fluid elements encounter one 
another and dissipate their energy. If this happens locally, then the ingredients 
for a classical a model of accretion are present. To the extent that many disk 
features are insensitive to, or independent of, the precise rms level of the 
disk fluctuations, these models can be of some practical utility. The classical 
formula relating the disk emissivity to radius R is the most important instance 
of this. 

In real protostellar disks, the level of ionization present is such that there 
are significant departures from the behavior of an ideal (perfectly conducting) 
MHD gas. In principle, the gas can become completely decoupled from the 
magnetic field, and go over to a hydrodynamic system. Less dramatically, the 
current-bearing electrons can acquire a velocity significantly different from 
the dominant neutrals, since the former need to maintain a current density 
and magnetic field even as charged species become rarer and rarer. When this 
happens, the field lines are no longer frozen into the motion of the bulk of 
the (neutral) disk gas, they are frozen into the electrons and the distinction 
becomes important. In the so-called Hall regime, the ions and neutrals move 
together distinctly from the electrons, and the ion motion relative to the 
electron-following field lines induces an additional electromotive force into 
the gas, beyond the self-induction responsible for simple field line freezing. 
Hall MHD is likely to be important in protostellar disks on scales of AU to 
lO's of AU scales. 

Dust grains are an ever present complication for MHD disk modeling. 
Typically, solids make up about 1 percent of the mass of interstellar gas. An 
interstellar population of small grains (radius ^ 10~^ cm) with a total mass 
fraction of a percent would present an enormous collecting area of would-be 
gas phase electrons. Putting charges on the grains effectively removes them 
because of the low mobility of the grains. Determining the disposition of the 
dust grains is thus a necessity for constructing MHD disk models. As the disk 
evolves, so too do the dust grains. They grow in size as they agglomerate, and 
they tend to settle towards the midplane, if they arc not stirred by turbulence 
or some other dynamical process. Larger dust grains are much less efficient in 
removing gas phase electrons than are smaller grains, where "larger" means 
growth in radius of an order of magnitude or more. 

Where does all of this leave us? Clearly, we are a long way from framing a 
picture of protostellar disk evolution at the level of, say, classical stellar evolu- 
tion. But for all of the gaping uncertainties, a useful zcroth order MHD-based 
picture of a protostellar disk in its T Tauri phase can be cobbled together: 

Inside of about 0.3 AU, a protostellar disk will be thermally ionized by 
direct exposure to the central source and self-consistcntly by the dissipation 
of MHD turbulence. We may expect vigorous MRI induced MHD turbulence 
in this zone. 

On scales of AU's, thermal ionization is no longer adequate to maintain 
the requiste levels of ionization to ensure MHD coupling, at least not near the 
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high density midplane. Depending upon the X-ray luminosity of the central 

star, the spectrum of dust grains and the abundance of gas phase potassium 
and sodium (electron donors), there could be good coupling at low density, 
higher elevation disk altitudes. The "dead zone" at lower altitudes need not be 
devoid of all transport; density waves from an adjacent active layer or large 
scale magnetic torques would each contribute their own stresses. It is even 
possible that some degree of coupling could be maintained down to the disk 
midplane, though this depends strongly on modeling assumptions. It seems 
likely, however, that on a scale an AU to tens of AU's , the level of MHD 
turbulence will be far less than in the disk's innermost regions. 

Finally, on scales larger than tens of AU's, the falling density lengthens 
recombination times and suggests a return to good MHD coupling. The typical 
size of a disk is many hundred AU's so that the dead zone is small when viewed 
globally, and may not have much of an impact in the overall accretion process 
(Terquem 2008). Perhaps, however, it is not coincidence that the MHD "quiet 
zone" coincides with the region of planet formation in the solar nebula and 
(more speculatively) in other protostellar disks. 
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Appendix A 

Begin with equation (fT6 



-XB=P7„+Pen (172) 

c 

where the right side of this equation is 

ninHln{<ynlWnl){vl - v) + Tig WTT-e (cr„e W„e) [(I'e ^ Vj) + (vj - v)] (173) 

In what foUows, we wiU often need an estimate of the ratio of the ion and 
electron coUision rates. FoUowing our discussion in section 3, we wiU assume 
that 

{(TneWne) \ e/^/n / 

where e < 1 is inserted because the electron-neutral cross section is geometri- 
cal, while the ion-ncutral collision cross section is larger than geometrical. 

In (|173p . the first {vi — v) term dominates the last by a factor of order 
(^/„/eTOe)^/^. Hence, 

J J r, nnie 

Pin = —XB - Pen = —XB -\ (cr„eWne)j (175) 

c c e 

With the help of equation (fTO)l . this implies 



J X B lerrie , , , 

VI -v = + J — -Z{vi-v,). (176) 

Cjppi V flnl 

A 

We have marked the term with an "A" for future reference. 
Recall equation ([^0]) : 

E+-[v + {ve- vi) + {vi -v)]xB+ ai^-^ [(ve - vi) + {vi - v)] = 0. 

c e 

(177) 

In substituting equation (|176p for vi — v in the above, we may always drop 
the "A" term, since it is small compared with — vi. 
Proceeding with the above substitution leads to 



V J X B 

E+-XB 



c erieC 



(JxB)xB J , , 



where (Jcond is defined in equation (|24p. The "B" term in the above equation 
may now clearly be dropped: it is of order Z(eTOeMn/)^^^/™n- This leads to 

^ V ^ JxB (J X B)xB J „ , 
E+-XB + ^ 5 '- = 0, (179) 

C eUeC C^JPPI (Jcond 



which is precisely equation (|23p in the text 
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Appendix B: 

To estimate the order of departures from charge neutrahty or the displacement 
currents, we wiU assume that the V operator is 1/^, where L is a charac- 
teristic length scale of the flow, and d/dt is ~ v/L, where w is a characteristic 
velocity (say, the largest of the neutral, ion, or electron velocities). 

For the electric field, we take E ~ vB/c, since we are only interested in 
problems where the inductive terms are comparable to, or in larger than, the 
resistive damping. Then 

V E - E/L - vB/Lc - AttvJ/c'^ - inen^v'^/c^ (180) 

Hence, the divergence of E is of order u^/c^ times the electron charge density 
(at most, since we assumed — We ~ w in the above). It may thus be ignored. 

For the displacement current, the demonstration is almost a matter of 
direct inspection. If we take the curl of equation ([^5]) and use equation ^7\ . 
then the on the right side of the first term is ^ B/L?, while the second, 
displacement, term is ~ v^Bjf? . It may thus be ignored. 



